Step 1: Define the problem

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
plt.style.use('seaborn')

import warnings
warnings.filterwarnings('ignore')
warnings.filterwarnings(action='ignore',category=DeprecationWarning)
warnings.filterwarnings(action='ignore',category=FutureWarning)

missing_values = ["?", ".", "", "_", "Na", "NULL", "null", "not", "Not", "NaN", "NA", "??", "nan", "inf"]
raw_data = pd.read_csv("Police_Department_Incident_Reports__2018_to_Present.csv", na_values=missing_values)
raw_data.set_index('Row ID', inplace=True)
raw_data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 406660 entries, 95300907041 to 95316305151
Data columns (total 25 columns):
 #   Column                   Non-Null Count   Dtype  
---  ------                   --------------   -----  
 0   Incident Datetime        406660 non-null  object 
 1   Incident Date            406660 non-null  object 
 2   Incident Time            406660 non-null  object 
 3   Incident Year            406660 non-null  int64  
 4   Incident Day of Week     406660 non-null  object 
 5   Report Datetime          406660 non-null  object 
 6   Incident ID              406660 non-null  int64  
 7   Incident Number          406660 non-null  int64  
 8   CAD Number               316015 non-null  float64
 9   Report Type Code         406660 non-null  object 
 10  Report Type Description  406660 non-null  object 
 11  Filed Online             83492 non-null   object 
 12  Incident Code            406660 non-null  int64  
 13  Incident Category        406376 non-null  object 
 14  Incident Subcategory     406376 non-null  object 
 15  Incident Description     406660 non-null  object 
 16  Resolution               406660 non-null  object 
 17  Intersection             385738 non-null  object 
 18  CNN                      385738 non-null  float64
 19  Police District          406660 non-null  object 
 20  Analysis Neighborhood    385652 non-null  object 
 21  Supervisor District      385738 non-null  float64
 22  Latitude                 385738 non-null  float64
 23  Longitude                385738 non-null  float64
 24  point                    385738 non-null  object 
dtypes: float64(5), int64(4), object(16)
memory usage: 80.7+ MB

raw_data.head()

	Incident Datetime	Incident Date	Incident Time	Incident Year	Incident Day of Week	Report Datetime	Incident ID	Incident Number	CAD Number	Report Type Code	...	Incident Description	Resolution	Intersection	CNN	Police District	Analysis Neighborhood	Supervisor District	Latitude	Longitude	point
Row ID
95300907041	15/08/2020 8:56	15/08/2020	8:56	2020	Saturday	15/08/2020 8:56	953009	200474239	NaN	VS	...	Vehicle, Recovered, Auto	Open or Active	NaN	NaN	Out of SF	NaN	NaN	NaN	NaN	NaN
95322706244	15/08/2020 9:40	15/08/2020	9:40	2020	Saturday	15/08/2020 18:21	953227	206121692	NaN	II	...	Theft, From Locked Vehicle, >$950	Open or Active	NaN	NaN	Park	NaN	NaN	NaN	NaN	NaN
64174871000	24/02/2018 22:00	24/02/2018	22:00	2018	Saturday	02/03/2018 10:13	641748	186051531	NaN	II	...	Lost Property	Open or Active	NaN	NaN	Park	NaN	NaN	NaN	NaN	NaN
95319604083	16/08/2020 3:13	16/08/2020	3:13	2020	Sunday	16/08/2020 3:14	953196	200491669	202290313.0	II	...	Firearm, Discharging in Grossly Negligent Manner	Open or Active	23RD ST \ ARKANSAS ST	23642000.0	Bayview	Potrero Hill	10.0	37.754827	-122.397729	(37.75482657770952, -122.39772873392515)
95326228100	16/08/2020 3:38	16/08/2020	3:38	2020	Sunday	16/08/2020 4:56	953262	200491738	202290404.0	II	...	Malicious Mischief, Breaking Windows	Open or Active	VALENCIA ST \ 15TH ST	24377000.0	Mission	Mission	9.0	37.766540	-122.422044	(37.76653957529556, -122.42204381448558)

5 rows × 25 columns

raw_data.describe()

	Incident Year	Incident ID	Incident Number	CAD Number	Incident Code	CNN	Supervisor District	Latitude	Longitude
count	406660.000000	406660.000000	4.066600e+05	3.160150e+05	406660.000000	3.857380e+05	385738.000000	385738.000000	385738.000000
mean	2018.876450	801205.474829	1.903172e+08	1.911781e+08	25104.691391	2.532790e+07	5.965443	37.769397	-122.423753
std	0.785939	103962.070832	8.805146e+06	1.901829e+07	25789.472280	3.085816e+06	2.786786	0.024048	0.026104
min	2018.000000	618687.000000	0.000000e+00	1.000000e+00	1000.000000	2.001300e+07	1.000000	37.707988	-122.511295
25%	2018.000000	711316.750000	1.808291e+08	1.824923e+08	6244.000000	2.397500e+07	3.000000	37.756352	-122.434036
50%	2019.000000	801363.500000	1.904171e+08	1.913912e+08	7046.000000	2.491800e+07	6.000000	37.775890	-122.417659
75%	2020.000000	891312.250000	2.000247e+08	2.002129e+08	61030.000000	2.642000e+07	8.000000	37.785829	-122.407337
max	2020.000000	981011.000000	9.811720e+08	1.000000e+09	75030.000000	5.412200e+07	11.000000	37.829991	-122.363743

Step 2: Data cleaning / preprocessing & Feature Engineering

Loại bỏ các biến thừa không cung cấp thông tin hữu ích nào cho việc phân tích

data = raw_data.drop([
    'Incident Date'       , 'Incident Time'      , 'Incident Year'          , 
    'Report Datetime'     , 'Report Type Code'   , 'Report Type Description',
    'Incident ID'         , 'Incident Number'    , 'CAD Number'             , 
    'Filed Online'        , 'Incident Code'      , 'Incident Subcategory'   , 
    'Incident Description', 'CNN'                , 'Supervisor District'    , 
    'point'
], axis=1)
data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 406660 entries, 95300907041 to 95316305151
Data columns (total 9 columns):
 #   Column                 Non-Null Count   Dtype  
---  ------                 --------------   -----  
 0   Incident Datetime      406660 non-null  object 
 1   Incident Day of Week   406660 non-null  object 
 2   Incident Category      406376 non-null  object 
 3   Resolution             406660 non-null  object 
 4   Intersection           385738 non-null  object 
 5   Police District        406660 non-null  object 
 6   Analysis Neighborhood  385652 non-null  object 
 7   Latitude               385738 non-null  float64
 8   Longitude              385738 non-null  float64
dtypes: float64(2), object(7)
memory usage: 31.0+ MB

Data Formating

data['Incident Datetime'] = pd.to_datetime(data['Incident Datetime'])
data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 406660 entries, 95300907041 to 95316305151
Data columns (total 9 columns):
 #   Column                 Non-Null Count   Dtype         
---  ------                 --------------   -----         
 0   Incident Datetime      406660 non-null  datetime64[ns]
 1   Incident Day of Week   406660 non-null  object        
 2   Incident Category      406376 non-null  object        
 3   Resolution             406660 non-null  object        
 4   Intersection           385738 non-null  object        
 5   Police District        406660 non-null  object        
 6   Analysis Neighborhood  385652 non-null  object        
 7   Latitude               385738 non-null  float64       
 8   Longitude              385738 non-null  float64       
dtypes: datetime64[ns](1), float64(2), object(6)
memory usage: 31.0+ MB

Loại bỏ các giá trị bị khuyết

data.isnull().sum().sort_values(ascending=False)

Analysis Neighborhood    21008
Longitude                20922
Latitude                 20922
Intersection             20922
Incident Category          284
Police District              0
Resolution                   0
Incident Day of Week         0
Incident Datetime            0
dtype: int64

data.dropna(axis=0, inplace=True)
data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 385370 entries, 95319604083 to 95316305151
Data columns (total 9 columns):
 #   Column                 Non-Null Count   Dtype         
---  ------                 --------------   -----         
 0   Incident Datetime      385370 non-null  datetime64[ns]
 1   Incident Day of Week   385370 non-null  object        
 2   Incident Category      385370 non-null  object        
 3   Resolution             385370 non-null  object        
 4   Intersection           385370 non-null  object        
 5   Police District        385370 non-null  object        
 6   Analysis Neighborhood  385370 non-null  object        
 7   Latitude               385370 non-null  float64       
 8   Longitude              385370 non-null  float64       
dtypes: datetime64[ns](1), float64(2), object(6)
memory usage: 29.4+ MB

Loại bỏ các giá trị trùng lặp

data.duplicated().sum()

24298

data.drop_duplicates(inplace=True)
data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 361072 entries, 95319604083 to 95316305151
Data columns (total 9 columns):
 #   Column                 Non-Null Count   Dtype         
---  ------                 --------------   -----         
 0   Incident Datetime      361072 non-null  datetime64[ns]
 1   Incident Day of Week   361072 non-null  object        
 2   Incident Category      361072 non-null  object        
 3   Resolution             361072 non-null  object        
 4   Intersection           361072 non-null  object        
 5   Police District        361072 non-null  object        
 6   Analysis Neighborhood  361072 non-null  object        
 7   Latitude               361072 non-null  float64       
 8   Longitude              361072 non-null  float64       
dtypes: datetime64[ns](1), float64(2), object(6)
memory usage: 27.5+ MB

Đặt tên lại cho các cột

data.rename(columns={
    'Incident Datetime': 'Datetime', 
    'Incident Day of Week': 'Weekday',
    'Incident Category': 'Category',
    'Police District': 'District',
    'Analysis Neighborhood': 'Neighborhood',
}, inplace=True)
data.info()

<class 'pandas.core.frame.DataFrame'>
Int64Index: 361072 entries, 95319604083 to 95316305151
Data columns (total 9 columns):
 #   Column        Non-Null Count   Dtype         
---  ------        --------------   -----         
 0   Datetime      361072 non-null  datetime64[ns]
 1   Weekday       361072 non-null  object        
 2   Category      361072 non-null  object        
 3   Resolution    361072 non-null  object        
 4   Intersection  361072 non-null  object        
 5   District      361072 non-null  object        
 6   Neighborhood  361072 non-null  object        
 7   Latitude      361072 non-null  float64       
 8   Longitude     361072 non-null  float64       
dtypes: datetime64[ns](1), float64(2), object(6)
memory usage: 27.5+ MB

Tạo một biến phân loại cho thời gian

def daypart(hour):
    if   hour in ['23','00','01','02']: return 'Midnight'
    elif hour in ['03','04','05','06']: return 'Early Morning'
    elif hour in ['07','08','09','10']: return 'Morning'
    elif hour in ['11','12','13','14']: return 'Noon'
    elif hour in ['15','16','17','18']: return 'Evening'
    return 'Night'

data['Daypart'] = data['Datetime'].dt.time.apply(lambda x: daypart(str(x).split(':')[0]))

data.head()

	Datetime	Weekday	Category	Resolution	Intersection	District	Neighborhood	Latitude	Longitude	Daypart
Row ID
95319604083	2020-08-16 03:13:00	Sunday	Assault	Open or Active	23RD ST \ ARKANSAS ST	Bayview	Potrero Hill	37.754827	-122.397729	Early Morning
95326228100	2020-08-16 03:38:00	Sunday	Malicious Mischief	Open or Active	VALENCIA ST \ 15TH ST	Mission	Mission	37.766540	-122.422044	Early Morning
95336264020	2020-08-16 13:40:00	Sunday	Non-Criminal	Open or Active	04TH ST \ MINNA ST	Southern	Financial District/South Beach	37.784044	-122.403712	Noon
95335012010	2020-08-16 16:18:00	Sunday	Weapons Offense	Cite or Arrest Adult	ORTEGA ST \ 48TH AVE	Taraval	Sunset/Parkside	37.751003	-122.507416	Evening
95300674000	2020-12-08 22:00:00	Wednesday	Missing Person	Open or Active	FILLMORE ST \ TURK ST	Northern	Western Addition	37.780496	-122.432140	Night

Phát hiện và loại bỏ outliers

plt.figure(figsize=(10, 5))

plt.subplot(1,2,1)
plt.boxplot(data['Latitude'])
plt.title('Latitude', fontsize=15)

plt.subplot(1,2,2)
plt.boxplot(data['Longitude'])
plt.title('Longitude', fontsize=15)

plt.show()

import scipy
for col in ['Latitude', 'Longitude']:
    prop = data[col]
    IQR = scipy.stats.iqr(prop)

    Q1 = np.percentile(prop, 25)
    Q3 = np.percentile(prop, 75)
    
    n_O_upper = data[prop > (Q3 + 1.5 * IQR)].shape[0]
    n_O_lower = data[prop < (Q1 - 1.5 * IQR)].shape[0]
    
    outliers_per = (n_O_upper + n_O_lower) / data.shape[0]
    print('Outliers Percentage of',  prop.name, ':', outliers_per)

Outliers Percentage of Latitude : 0.007117693977932379
Outliers Percentage of Longitude : 0.06820800283599947

=> Có thể không cần loại bỏ Outliers vì số lượng không đáng kể

Exploratory Data Analysis (EDA)

plt.figure(figsize=(10, 7))

top10_category = data['Category'].value_counts()[:10]
ax = sns.countplot(
    y = 'Category', 
    data = data, 
    order = top10_category.index
)

for rect in ax.patches:
    ax.text(
        rect.get_width() - 8000, 
        rect.get_y() + rect.get_height() / 2,
        rect.get_width(),
        color = 'white',
        weight = 'bold',
        fontsize = 11
    )
    
plt.title('Top 10 loại tội phạm ở San Francisco', fontsize=16)
plt.show()

=> Tội phạm thuộc loại Larceny Theft, cao hơn đáng kể so với bất kỳ loại tội phạm nào khác

from matplotlib.animation import FuncAnimation
from IPython.display import HTML
from datetime import datetime

df = data.copy()
df['Year'] = df['Datetime'].dt.year

year_counts = df.groupby(['Year', 'Category']).count().reset_index().pivot(
    index = 'Year', 
    columns = 'Category', 
    values = 'Daypart'
)
year_counts.fillna(0, inplace=True)

def style_axes(ax):
    ax.tick_params(labelsize=5, length=0)
    ax.grid(True, axis='x', color='white')
    ax.set_axisbelow(True)
    [spine.set_visible(False) for spine in ax.spines.values()]
    
def prepare_data(df, steps=20):
    df = df.reset_index()
    df.index = df.index * steps
    
    df_expanded = df.reindex(range(df.index[-1] + 1))
    df_expanded['Year'] = df_expanded['Year'].fillna(method='ffill')
    
    df_expanded = df_expanded.set_index('Year')
    df_rank_expanded = df_expanded.rank(axis=1, method='first')
    
    df_expanded = df_expanded.interpolate()
    df_rank_expanded = df_rank_expanded.interpolate()
    
    return df_expanded, df_rank_expanded

def init():
    ax.clear()
    style_axes(ax)

def update(i):
    ax.clear()
    for bar in ax.containers: bar.remove()
    y = df_rank_expanded.iloc[i]
    
    ax.barh(
        y = y, 
        width = df_expanded.iloc[i], 
        color = plt.cm.Dark2(range(6)), 
        tick_label = df_expanded.columns
    )
    
    for rect in ax.patches:
        ax.text(
            rect.get_width() + 500, 
            rect.get_y() + rect.get_height() / 4,
            int(rect.get_width()),
            color = 'blue',
            fontsize = 5.5
        )
    
    ax.set_ylim(min(y) - 1, max(y) + 1)
    ax.set_title(
        f'Tội phạm ở San Francisco - Năm {int(df_expanded.index[i])}', 
        fontsize=10
    )

fig = plt.Figure(figsize=(6.8, 6), dpi=144)
ax = fig.add_subplot()

df_expanded, df_rank_expanded = prepare_data(year_counts)
animation = FuncAnimation(
    fig = fig, 
    func = update, 
    init_func = init, 
    frames = len(df_expanded), 
    interval = 100, 
    repeat = False
)
HTML(animation.to_jshtml())

Once Loop Reflect

Phân tích theo Resolution

color = plt.cm.winter(np.linspace(0, 10, 20))
data['Resolution'].value_counts().plot.bar(color=color)

plt.xticks(rotation=0, fontsize=11)
plt.title('Các giải pháp cho tội phạm',fontsize=16)
plt.show()

=> Hầu hết các sự cố đều đang được tiến hành xử lý

most_commons = data[data['Category'].isin(top10_category.index)]
violent = most_commons.copy()
violent['Arrest'] = np.where(violent['Resolution'].isin(['Cite or Arrest Adult', 'Unfounded']), 0, 1)

arrest_counts = violent['Category'][violent.Arrest == 1].value_counts()[:10]
total_counts = violent['Category'].value_counts()[:10]
arrest_counts = arrest_counts / (total_counts).sort_index()
total_counts = total_counts / (total_counts).sort_index()

total_counts.plot.barh(color='crimson', label= 'Unsolved')
arrest_counts.plot.barh(color='mediumseagreen', label='Solved')

plt.legend(bbox_to_anchor=(1.05, 1))
plt.show()

plt.figure(figsize = (12, 10))

arrest = data[data['Resolution'].isin(['Cite or Arrest Adult', 'Cite or Arrest Juvenile'])]
a = arrest['Category'].value_counts().reset_index()
b = data['Category'].value_counts().reset_index()

a = a.merge(b, how='inner', on='index')
a.columns = ['Category', 'Arrests', 'Cases']

a['Arrests Percent'] = round(100 * (a['Arrests'] / a['Cases']), 2)
a.sort_values('Arrests Percent', ascending = False, inplace=True)

sns.barplot(x='Arrests Percent', y='Category', data=a)
plt.xlabel('')
plt.ylabel('')

plt.title('Điều gì sẽ khiến ta bị bắt ở San Francisco', fontsize=16)
plt.show()

=> Ma túy, vi phạm giao thông và lệnh bắt giữ là 3 lý do hàng đầu cho các vụ bắt giữ

Phân tích theo Quận

plt.figure(figsize=(10, 7))
color = plt.cm.spring(np.linspace(0, 1, 12))
ax = data['District'].value_counts().plot.bar(color=color)

for p in ax.patches:
    height = p.get_height()
    ax.text(
        x = p.get_x() + p.get_width() / 2, 
        y = height + 1000, 
        s = height, 
        ha = 'center'
    )
    
plt.xticks(rotation=0, fontsize=11)
plt.title('Quận có nhiều tội phạm nhất', fontsize=16)
plt.show()

=> Quận Central là nơi ghi nhận nhiều tội phạm nhất và quận Park là nơi ít tội phạm nhất ở San Francisco

df = pd.crosstab(data['Category'], data['District'])
color = plt.cm.Greys(np.linspace(0, 1, 10))
df.div(df.sum(1).astype(float), axis=0).plot.bar(stacked=True, color=color, figsize=(15, 7))

plt.legend(loc='upper left', bbox_to_anchor=(1, 0, 0.5, 1), fontsize=12)
plt.xlabel('')

plt.title('Quận vs Loại tội phạm', fontsize=16)
plt.show()

plt.figure(figsize = (12, 5))

a = arrest['District'].value_counts().reset_index()
b = data['District'].value_counts().reset_index()

a = a.merge(b,how = 'inner', on = 'index')
a.columns = ['District', 'Arrests','Cases']

a['Arrests Percent'] = round(100 * a['Arrests'] / a['Cases'], 2)
a.sort_values('Arrests Percent', ascending = False, inplace = True)

sns.barplot(
    x = 'District', 
    y = 'Arrests Percent',
    data = a, 
    palette = ['red'] + ['grey'] * 8 + ['blue'] + ['grey']
)

for i in range(10):
    plt.text(
        x = i, 
        y = a['Arrests Percent'] .iloc[i], 
        s = f"{a['Arrests Percent'].iloc[i]}%",
        horizontalalignment = 'center', 
        verticalalignment = 'bottom',
        weight = 'bold'
    )

plt.xticks(fontsize=11)    
plt.xlabel('')
plt.ylabel('')

plt.title('Nơi mà cảnh sát nghiêm khắc nhất - Phần trăm bị bắt', fontsize=16)
plt.show()

=> Mặc dù có ít tội phạm hơn được ghi nhận ở quận Tenderlion nhưng tỷ lệ bắt giữ ở đây khá cao 36,13%

Phân tích theo Thời gian

del df
df = data.copy()
df = df.groupby(df['Datetime'].dt.date).count().iloc[:, 0]

color_palette = sns.color_palette()
sns.kdeplot(data=df, shade=True)

plt.axvline(
    x = df.median(), 
    ymax = 0.95, 
    linestyle = '--', 
    color = color_palette[1]
)

plt.annotate(
    f'Median: {df.median()}',
    xy = (df.median(), 0.004),
    xytext = (200, 0.005),
    arrowprops = dict(
        arrowstyle = '->', 
        color = color_palette[1], 
        shrinkB = 10
    )
)

plt.xlabel('Incidents')
plt.ylabel('Density')

plt.title('Phân phối số lượng sự cố mỗi ngày', fontdict={'fontsize': 16})
plt.show()

cat_per_week_common = pd.crosstab(most_commons['Category'], most_commons['Weekday']) 
cat_per_week_common = cat_per_week_common.div(cat_per_week_common.sum(axis=1), axis=0)
cat_per_week_common = cat_per_week_common[['Monday', 'Tuesday', 'Wednesday', 'Thursday', 'Friday', 'Saturday','Sunday']]

ax = sns.heatmap(cat_per_week_common, cmap="BuPu", linewidths=0.5)
plt.xticks(fontsize=11, rotation=45, ha='right')
plt.yticks(fontsize=11)

plt.xlabel('')
plt.ylabel('')
plt.show()

data['Weekday'].value_counts().plot.pie(
    explode = [0.1, 0, 0, 0, 0, 0, 0], 
    textprops = {'color': 'w', 'fontsize': 15, 'weight': 'bold'},
    autopct = "%.1f%%"
)

plt.legend(loc='center left', bbox_to_anchor=(1, 0, 0.5, 1), fontsize=12)
plt.axis('off')

plt.title('Số lượng tội phạm theo thứ', fontsize=16)
plt.show()

=> Thứ 6 là ngày mà tội phạm được ghi nhận nhiều nhất, tiếp theo là Thứ 4. Chủ nhật là ngày ít tội phạm nhất

del df
df = data.copy()
df['Date'] = df['Datetime'].dt.date
df['Hour'] = df['Datetime'].dt.hour

df = df.groupby(['Hour', 'Date', 'Category'], as_index=False).count().iloc[:, :4]
df.rename(columns={'Datetime': 'Incidents'}, inplace=True)

df = df.groupby(['Hour', 'Category'], as_index=False).mean()
df = df[df['Category'].isin(top10_category[4:].index)]

fig, ax = plt.subplots(figsize=(12, 5))
ax = sns.lineplot(x='Hour', y='Incidents', data=df, hue='Category')
ax.legend(loc='upper center', bbox_to_anchor=(0.5, 1.15), ncol=6)
fig.tight_layout(rect=[0, 0.03, 1, 0.95])

plt.suptitle('Số sự cố trung bình mỗi giờ')
plt.show()

plt.figure(figsize=(15, 7))
color = plt.cm.twilight(np.linspace(0, 4, 100))
data['Datetime'].dt.time.value_counts().head(24).sort_index().plot.bar(color=color)

plt.xticks(rotation=45, fontsize=11)
plt.title('Phân bố tội phạm trong ngày', fontsize=16)
plt.show()

df = data['Daypart'].value_counts().reset_index()
df_normalize = data['Daypart'].value_counts(normalize=True)
sns.barplot(x=df['index'], y='Daypart', data=df)

for i in range(6):
     plt.text(
         x = i, 
         y = 5000, 
         s = f'{round(100 * df_normalize[i], 2)}%', 
         horizontalalignment = 'center', 
         color = 'white',
         weight = 'bold'
     )
        
plt.title('Số lượng tội phạm theo các buổi trong ngày', fontsize=16)
plt.show()

=> Hãy cẩn thận vào buổi tối. Gần 1/4 số vụ tội phạm xảy ra vào buổi tối

plt.figure(figsize = (15, 7))

pivot_table = pd.pivot_table(
    columns = data['Daypart'] ,
    index = 'Weekday',
    values = 'Daypart' , 
    aggfunc = 'count', 
    data = data
)

pivot_table = pivot_table.reindex(
    index = ['Monday', 'Tuesday', 'Wednesday', 'Thursday', 'Friday', 'Saturday', 'Sunday'], 
    columns = ['Early Morning', 'Morning', 'Noon', 'Evening', 'Night', 'Midnight']
)

sns.heatmap(pivot_table, cmap='Reds', annot=True, fmt='d')
plt.xlabel('')
plt.ylabel('')
plt.show()

=> Buổi tối thường là lúc có nhiều vụ phạm tội, tối Thứ 6 đặc biệt tồi tệ và tối Chủ nhật tương đối tốt hơn

sns.countplot(data['Datetime'].dt.month, palette='autumn')
plt.xticks(fontsize=11)
plt.xlabel('')

plt.title('Tội phạm theo từng tháng', fontsize=16)
plt.show()

=> Tháng 1 là tháng có nhiều tội phạm nhất

plt.figure(figsize = (15, 7))

df = data[['Datetime']]
df['Day'] = df['Datetime'].dt.day
df['Month'] = df['Datetime'].dt.month

pivot_table = pd.pivot_table(
    columns = df.Day, 
    values = 'Day', 
    index = 'Month', 
    aggfunc = 'count', 
    data = df
)
sns.heatmap(pivot_table, cmap='Reds')

plt.xticks(fontsize=12)
plt.yticks(
    np.arange(0.5, 12.5),
    ['Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sept', 'Oct', 'Nov', 'Dec'],
    fontsize = 12
)

plt.xlabel('')
plt.ylabel('')
plt.show()

=> Ngày 1/1 dường như là ngày có số tội phạm được ghi nhận cao nhất và ngày 25/12 là ngày có số tội phạm được ghi nhận ít nhất. Ngoài ra có thể thấy ngày 12 trong năm là ngày khá an toàn so với phần lớn các ngày còn lại

Phân tích Time Series

first_date = data['Datetime'].min()
last_date = data['Datetime'].max()
top10_category

Larceny Theft          107521
Other Miscellaneous     26798
Non-Criminal            23929
Malicious Mischief      23807
Assault                 22031
Burglary                19565
Motor Vehicle Theft     16510
Warrant                 12704
Fraud                   10935
Lost Property           10584
Name: Category, dtype: int64

print('Ngày đầu tiên:', first_date)
print('Ngày cuối cùng:', last_date)
print(f'Số ngày tổng cộng: {(last_date - first_date).days + 1}')

Ngày đầu tiên: 2018-01-01 00:00:00
Ngày cuối cùng: 2020-12-11 23:58:00
Số ngày tổng cộng: 1076

time_span = pd.date_range(first_date, last_date)
def to_time_series(df):
    date = df['Datetime'].value_counts().index
    ts = pd.DataFrame(
        data = df['Datetime'].value_counts().values[np.argsort(date)], 
        index = sorted(date), 
        columns = ['count']
    )
    return ts.reindex(time_span, fill_value=0) 

all_counts = to_time_series(data)
all_counts[all_counts['count'] == 0]

	count
2020-01-12	0
2020-02-12	0
2020-03-12	0
2020-04-12	0
2020-05-12	0
2020-06-12	0
2020-07-12	0
2020-08-12	0
2020-09-12	0
2020-10-12	0
2020-11-12	0
2020-11-19	0
2020-11-23	0
2020-11-24	0
2020-11-25	0
2020-11-26	0
2020-11-27	0
2020-11-28	0
2020-11-29	0
2020-11-30	0

=> Vậy những ngày này liệu có phải là những ngày bình yên hay đây có lẽ chỉ là sơ suất nào đó của cơ quan cảnh sát. Từ biểu đồ trước đó cùng với phân tích này, ta có thể nhận thấy rõ điều có gì đó ở ngày 12 khiến cho tội phạm vào ngày này khá thấp

plt.figure(figsize=(16, 7))
all_counts['count'].plot()

plt.xlabel('Thới gian')
plt.ylabel('Số lượng tội phạm')
plt.show()

=> Có vẻ như tỷ lệ tội phạm không thay đổi nhiều trong những năm qua! Liệu nó có đúng cho các loại riêng lẻ? Ta sẽ thử với Larceny Theft

df_theft = data[data['Category'] == 'Larceny Theft']
theft_counts = to_time_series(df_theft)
theft_counts[theft_counts['count'] == 0]

	count
2018-02-11	0
2018-02-22	0
2018-03-01	0
2018-04-12	0
2018-05-02	0
...	...
2020-11-27	0
2020-11-28	0
2020-11-29	0
2020-11-30	0
2020-12-02	0

84 rows × 1 columns

plt.figure(figsize=(16, 7))
theft_counts['count'].plot()

plt.xlabel('Thới gian')
plt.ylabel('Số lượng tội phạm Larceny Theft')
plt.show()

=> Nó thực sự không thay đổi gì nhiều

Phân tích theo Mô tả

data_with_description = raw_data[raw_data.index.isin(data.index)]
data_with_description.rename(columns={
    'Incident Datetime': 'Datetime', 
    'Incident Day of Week': 'Weekday',
    'Incident Category': 'Category',
    'Incident Description': 'Description',
    'Police District': 'District',
    'Analysis Neighborhood': 'Neighborhood',
}, inplace=True)

from wordcloud import WordCloud
description = data_with_description['Description']
wc = WordCloud(background_color='#e9eaf1', width=1000, height=500).generate(str(description))

plt.imshow(wc)
plt.axis('off')

plt.title('Mô tả cho các tội', fontsize=20)
plt.show()

import plotly.express as px
discat = data_with_description.groupby(['District', 'Category'])['Description'].count().reset_index()
fig = px.sunburst(
    discat.rename(columns={'Description': 'Count'}), 
    path = ['District', 'Category'], 
    values = 'Count', 
    color = 'Count'
)
fig.show()

Phân tích theo Phân bố không gian (Spatial Distribution)

plt.figure(figsize=(10, 7))
ax = sns.countplot(
    y = 'Intersection', 
    data = data, 
    order = data['Intersection'].value_counts()[:15].index
)

for rect in ax.patches:
    ax.text(
        rect.get_width() - 150, 
        rect.get_y() + rect.get_height() / 2,
        rect.get_width(),
        color = 'white',
        weight = 'bold', 
        fontsize = 11
    )
    
plt.title('Top 15 vùng phạm tội', fontsize=16)
plt.show()

del df
df = data.copy()

mean_lat = np.mean(df.Latitude)
mean_lon = np.mean(df.Longitude)

sw = df[['Latitude', 'Longitude']].min().values.tolist()
ne = df[['Latitude', 'Longitude']].max().values.tolist()

df['Coordinates'] = df['Longitude'].astype(str) + ', ' + df['Latitude'].astype(str)
df.Coordinates.value_counts()[:5]

-122.40733700000001, 37.78456014    2484
-122.4036355, 37.77516081           2278
-122.4080362, 37.78640961           2157
-122.419669, 37.76505134            1666
-122.41259529999999, 37.78393258    1488
Name: Coordinates, dtype: int64

=> Cặp tọa độ phổ biến nhất (-122.40733700000001, 37.78456014) là vị trí của Market St

plt.title('SF Crime Distribution')
plt.xlabel('Longitude')
plt.ylabel('Latitude')
plt.scatter(data.Longitude, data.Latitude, s=0.5, c='r')
plt.show()

=> Ta có thể thấy có nhiều tội phạm hơn ở phần đông bắc của San Francisco. Câu hỏi đặt ra rằng liệu điều này có luôn như vậy trong suốt các năm qua không ?

Ta sẽ xem liệu phân phối có bất kỳ thay đổi đáng kể nào không:

• Để giảm thời gian tính toán, ở đây sẽ chỉ vẽ dữ liệu của tháng 1 (tháng có nhiều tội phạm nhất) hàng năm.
• Để hiểu rõ hơn về sự phân bố không gian, lần này biểu đồ sẽ được phủ các đường bao lên các ô phân tán.

plt.figure(figsize = (16, 5))
for index, year in enumerate(range(2018, 2021)): 
    df_jan = data[(data['Datetime'] < f'{year}-02') & (data['Datetime'] > f'{year - 1}-12-31')]
    
    plt.subplot(1, 3, index + 1)
    plt.scatter(df_jan.Longitude, df_jan.Latitude, s=0.5, c='r')
    sns.kdeplot(df_jan.Longitude, df_jan.Latitude)
    
    plt.xticks(rotation=45)
    plt.xlabel('Longitude')
    plt.ylabel('Latitude')
    
    plt.title(f'Phân bố tội phạm vào tháng 1 năm {year}')
plt.show()

=> Có vẻ như phần đông bắc của San Francisco luôn là khu vực nguy hiểm nhất kể từ năm 2018. Ngoài ra, các đường nét cho thấy sự phân bố tội phạm xung quanh khu vực này đã mở rộng 1 chút về phía đông kể từ năm 2020

Trực quan hóa không gian địa lý

import geopandas as gpd
from shapely.geometry import Point

gdf = data.copy()
gdf['Coordinates'] = list(zip(gdf.Longitude, gdf.Latitude))
gdf.Coordinates = gdf.Coordinates.apply(Point)

gdf = gpd.GeoDataFrame(gdf, geometry='Coordinates', crs={'init': 'epsg:4326'})
world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))

ax = world.plot(color='white', edgecolor='black')
gdf.plot(ax=ax, color='red')
plt.show()

district_labels = {
    'Southern': '37.774432, -122.401121',
    'Bayview': '37.734332, -122.389920',
    'Mission': '37.756478, -122.423663',
    'Northern': '37.787740, -122.430300',
    'Tenderloin': '37.781980, -122.412981',
    'Central': '37.796200, -122.409293',
    'Park': '37.765352, -122.449282',
    'Richmond': '37.776204, -122.483285',
    'Ingleside': '37.726817, -122.437207',
    'Taraval': '37.737775, -122.484375'
}

table = data['District'].value_counts().reindex(district_labels.keys())
table = table.reset_index().rename(
    {'index': 'District', 'District': 'Count'}, 
    axis='columns'
)

table['District'] = table['District'].str.upper()
table

	District	Count
0	SOUTHERN	44566
1	BAYVIEW	31776
2	MISSION	48478
3	NORTHERN	50401
4	TENDERLOIN	34886
5	CENTRAL	54642
6	PARK	17457
7	RICHMOND	21149
8	INGLESIDE	28019
9	TARAVAL	25575

import urllib.request
import shutil
import zipfile

url = 'https://data.sfgov.org/api/geospatial/wkhw-cjsf?method=export&format=Shapefile'
with urllib.request.urlopen(url) as response, open('pd_data.zip', 'wb') as out_file:
    shutil.copyfileobj(response, out_file)
    
with zipfile.ZipFile('pd_data.zip', 'r') as zip_ref:
    zip_ref.extractall('pd_data')

import os
import re
    
for filename in os.listdir('./pd_data/'):
    if re.match(".+\.shp", filename):
        districts = gpd.read_file(f'./pd_data/{filename}')
        break

districts.crs={'init': 'epsg:3857'}
districts = districts.merge(
    table.set_index(['District']),
    how = 'inner',
    left_on = 'district',
    right_index = True,
    suffixes = ('_x', '_y')
)
districts

	company	district	shape_area	shape_le_1	shape_leng	geometry	Count
0	B	SOUTHERN	9.134414e+07	100231.353916	87550.275142	MULTIPOLYGON (((-122.39186 37.79425, -122.3917...	44566
1	C	BAYVIEW	2.013846e+08	144143.480351	163013.798332	POLYGON ((-122.38098 37.76480, -122.38103 37.7...	31776
2	D	MISSION	8.062384e+07	40518.834235	40152.783389	POLYGON ((-122.40954 37.76932, -122.40862 37.7...	48478
3	E	NORTHERN	8.278169e+07	50608.310321	56493.858208	POLYGON ((-122.43379 37.80793, -122.43375 37.8...	50401
4	J	TENDERLOIN	1.107215e+07	18796.784185	12424.268969	POLYGON ((-122.40217 37.78626, -122.41718 37.7...	34886
5	A	CENTRAL	5.595027e+07	67686.522865	64025.129073	POLYGON ((-122.42612 37.80684, -122.42612 37.8...	54642
6	F	PARK	8.487896e+07	50328.913294	46307.776968	POLYGON ((-122.43956 37.78314, -122.43832 37.7...	17457
7	G	RICHMOND	1.379640e+08	75188.628361	69991.465355	POLYGON ((-122.44127 37.79149, -122.44060 37.7...	21149
8	H	INGLESIDE	1.935805e+08	74474.181164	74737.936295	POLYGON ((-122.40450 37.74858, -122.40407 37.7...	28019
9	I	TARAVAL	2.846767e+08	73470.424000	75350.217521	POLYGON ((-122.49842 37.70810, -122.49842 37.7...	25575

import contextily as ctx
districts['incidents_per_day'] = districts.Count / data.groupby('Datetime').count().shape[0]

fig, ax = plt.subplots(figsize=(12, 7))
districts.plot(
    column = 'incidents_per_day',
    cmap = 'Reds',
    alpha = 0.6,
    edgecolor = 'r',
    linestyle = '-',
    linewidth = 1,
    legend = True,
    ax = ax)

for index in districts.index:
    name = districts.loc[index].district
    geometry = districts.loc[index].geometry
    
    plt.annotate(
        name,
        (geometry.centroid.x, geometry.centroid.y),
        color = '#353535',
        fontsize = 'large',
        fontweight = 'heavy',
        horizontalalignment = 'center'
    )
plt.show()

import geoplot as gplt
crimes = data['Category'].unique().tolist()

sf_land = districts.unary_union
sf_land = gpd.GeoDataFrame(gpd.GeoSeries(sf_land), crs={'init':'epsg:4326'})
sf_land = sf_land.rename(columns={0: 'geometry'}).set_geometry('geometry')

fig, ax = plt.subplots(2, 3, sharex=True, sharey=True, figsize=(15, 7))
for i, crime in enumerate(np.random.choice(crimes, size=6, replace=False)):
    ax = fig.add_subplot(2, 3, i + 1)
    gplt.kdeplot(
        gdf.loc[gdf['Category'] == crime],
        shade = True,
        shade_lowest = False,
        clip = sf_land.geometry,
        cmap = 'Reds',
        ax = ax)
    gplt.polyplot(sf_land, ax=ax)
    ax.set_title(crime) 

fig.tight_layout(rect=[0, 0, 0, 0])
ax.set_axis_off()

plt.title('Geographic Density of Different Crimes')
plt.show()

import folium
gjson = r'https://cocl.us/sanfran_geojson'

figure = folium.Figure(width=900, height=600)
sf_map = folium.Map(location=[37.77, -122.42], zoom_start=12)
sf_map.choropleth(
    geo_data = gjson,
    data = table,
    columns = ['District', 'Count'],
    key_on = 'feature.properties.DISTRICT',
    fill_color = 'YlOrRd', 
    fill_opacity = 0.7, 
    line_opacity = 0.2,
    legend_name = 'Crime Rate in San Francisco'
)

figure.add_child(sf_map)
figure

%%html
<iframe 
    src="https://data.sfgov.org/dataset/Map-of-Police-Department-Incident-Reports-2018-to-/jq29-s5wp/embed?width=950&height=600" 
    width="950" 
    height="600" 
    style="border:0; padding: 0; margin: 0;"
></iframe>

Step 3 : Model selection & Trainning

Step 4: Evaluate & test the model results & performance

Step 5: Deploy the model